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We develop a model in two dimensions to characterise the growth rate of a tracer 

^^ gradient mixed by a statistically homogeneous flow with rapid temporal variations. 

f--N The model is based on the orientation dynamics of the passive-tracer gradient with 

CN respect to the straining (compressive) direction of the flow, and involves reducing 

^^ the dynamics to a set of stochastic differential equations. The statistical properties 

u of the system emerge from solving the associated Fokker-Planck equation. In a 

^ certain limiting case, and within the model framework, there is a rigorous proof 

'nT that the tracer gradient aligns with the compressive direction. This limit involves 

I— I decorrelated flows whose mean vorticity is zero. Using numerical simulations, we 

^ll assess the extent to which our model applies to real mixing protocols, and map the 

'^ stochastic parameters on to flow parameters. 

I 

j5 PACS numbers: 47.27.wj, 05.10.Gg 
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>. I. INTRODUCTION 

'— ' When a passive tracer is stirred by a flow that varies rapidly in space and time, stochastic 

,__! models are often used to describe the mixing pj. Based on the Fokker-Planck (FP) equation, 

^ we develop a model for the probability distribution function (PDF) of the growth rate of 

the tracer gradient in a generic, two-dimensional, rapidly-varying stirring flow. In this 
l/~j Introduction, we place our model in the context of existing literature. 

O Kraichnan [2j modelled homogeneous isotropic turbulence as Gaussian fluctuations. In 

•n this description, the PDF for the largest finite-time Lyapunov exponent (FTLE) can be found 

using the Central Limit Theorem. The same model has been used to describe the so-called 
strange eigenmode [3J in passive-tracer decay [1] . Balkovsky and Fouxon [5] generalised the 
Kraichnan model to account for temporally-correlated flows. When the correlation time 
tends to zero, their analysis yields the PDF of all the FTLEs (in D dimensions) [5* iQ • For a 
review of these techniques, see [T]. In the present work, we focus instead on the distribution 
of the growth rate of the tracer gradient, a distinct but related quantity. 

Lapeyre et al. [7] derived a pair of equations to describe the growth rate of the gradient 
of a passively-advected tracer. The theory is based on the orientation of the tracer gradient 
relative to the straining (compressive) direction of the flow: if the tracer gradient aligns with 
the compressive direction, it blows up. However, the tracer-gradient orientation is instead 
flxed by a balance between strain and 'effective' vorticity (the sum of vorticity and twice 
the rotation rate of the strain eigenbasis). Thus, blowup is prevented in rotation-dominated 
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regions. The inclusion of the rotation rate generahses the description of orientation dynam- 
ics imphed by the Okubu- Weiss criterion [8]. That said, the Lapeyre criterion rehes on 
the adiabatic approximation of the differential equations: the effective vorticity and strain 
eigenvalue are assumed to vary slowly along trajectories. Although this description is valid 
in certain cases (e.g. [9]), it may not hold when the stirring flow is governed by some kind 
of externally-imposed, rapidly-varying forcing: this will be our focus herein. 

Gonzalez and Paranthoen [10] studied two mixing protocols wherein the adiabatic as- 
sumption breaks down. They show numerically that the tracer gradient points a direction 
fixed by the average value of the ratio r of effective vorticity to strain. Garcia and cowork- 
ers [11] introduced a stochastic model of such scenarios by regarding r as the solution to a 
stochastic differential equation (SDE). They found good agreement with an experiment they 
performed. However, they neither examined the distribution of growth rates of the tracer 
gradient elicited by this model (but see [12]), nor analysed the associated FP equation. In 
the present work, we address these issues and provide a more detailed model of a generic 



flow (Sec. III). We make the connection with the FP equation [13], an approach that facili- 
tates fast, accurate numerical calculations, analytical results in a certain limiting case, and 
highlights the striking analogies between the current work and a variety of other physical 
systems [13]. Furthermore, our approach involves a comparison between theoretical models 



and simulations of the random-phase sine flow (Sec. IV), and forced, two-dimensional tur- 
bulence (Sec. Ivl), an exercise that strengthens the case for model developed herein. First, 
we review the orientation dynamics derived elsewhere [71 [H] . 



II. THE ORIENTATION DYNAMICS 

The description of the orientation dynamics starts with the vector field S = {—6y,6x), 
where 9 is the passively-advected tracer 

9t + u-Ve = kAO, V ■u = 0. (1) 

In this section, we take the molecular diffusivity k, to be zero; we study the implications 
of finite K-values at the end of the derivation. When k = 0, the vector 13 satisfies the 
material-line equation 

Bt + uVB = BVu. (2) 

This result can be obtained by brute- force calculation based on Eq. ([T|, or by other 
means [I5] ^. 

>^W.B.)(^t)(|). (3) 

The rate-of-strain matrix S appears in Eq. ([3]): s = u^ and d = (uy + v^) /2. We identify 
/3 — (7r/2) as the phase of the vector B, where tan/3 = 6y/6x, and obtain, by means similar 
to before, an equation for /3: 



fH-^(«.-«.)(rM(|i. (4) 



^ Indeed, one may regard B as a complex-valued two-form B = Bi + \B2 and operate on both sides of the 
two-form identity Bdz A dz = 2idz A A9 with the material derivative. 



where u = Vx — Uy is the vorticity. To complete our analysis, we re-write Eq. ([s]) in terms of 
the eigenvalues of S: 

^|B|2 = -2AsinC|Bp. (5) 

Here, A = sign(d)\/s^~+~(P is the unsigned eigenvalue of 5*, with associated eigenvector X(+). 



The angle between the x-axis and ^(+) is </?, and ( = 2{(3 — if + ^tt). Hence, from Eq. (|4 

^ = -2AcosC + w + 2^. (6) 

dt dt ^ ^ 

Finally, we identify the growth rate of the gradient: 

A:=-2AsinC. (7) 

To understand this description, imagine a situation where A > 0. Then the growth rate is 
maximised when ( = — 7r/2, or P = ip — tt/2, that is, when the tracer gradient aligns with 
the compressive direction of the flow. In this paper, our concern is the extent to which this 
alignment condition is fulfilled. 

The basic equations (tsl) and ^ rely only on the assumptions k = and incompressibility. 
If K is finite but small, then the main effect of diffusion will be to reduce the growth rate 
of the gradient, through the presence of a negative term in Eq. ([s]). We do not expect it to 
interfere much with the orientation dynamics [3, [H]. Hence, we expect similar results for 
the PDF of — 2A sin (, with or without a small amount of diffusion. This is confirmed in our 
numerical simulations in Sec. |Vl 



III. A STOCHASTIC MODEL 

For notational reasons, we write X := ( [13] and re- write Eq. dol): 

dX 
\— = (-AocosX + w)-Y (t) cosX + Zo (t) , (8) 

where A = Xq + Y (t) and (w/2) + = w + Zo(t) represent decompositions into mean 
components and fluctuations. We assume that the vector X(+) spends an equal amount of 
time representing the compressive and expansive directions, hence Aq = 0. We model the 
fluctuations as Ornstein-Uhlenbeck (OU) processes, with mean zero, common decay time r, 
and strengths Dy and Dz respectively (the fluctuations reduce to Wiener processes in the 
limit where the decay time tends to zero). We anticipate this model to be appropriate in 
particular when the flow is governed by some random, externally-prescribed stirring mech- 
anism (for examples, see [IE]). For generality, we assume that the noise terms in Eq. ([s]) 
satisfy the following correlation relations: 

(r(t)F(t')) = :^e-l*-*'l/^ 

r 

(Zo(t)Zo(t')) = ^e-l*-*'l/^ 
r 

(r(t)Zo(t')) = c^^^5^e-l*-'l/^ 



where < c < 1 is a correlation coefficient, and where r~^e~l*l/'^ converges in the sense 
of distributions to 25 (s) as r — )■ 0. Finally, we make the assumption that the underlying, 
noise-generating flow is homogeneous in space in a statistical sense: 

(F it- x^) Y it'; x',)) = t-'Vy {x, - x',) e^^'-''^'^ &c. 

Hence, 

(F (t; xo) Y it'- xo)) = t'^Vy (0) e-l*-*'l/- := r-^Dye-l^^^'l/^ &c., 

and the noise strength is independent of the initial position of the Lagrangian particle. 

Equation rtsl) is Markovian with uncorrelated noise terms if viewed in an augmented state 
space: 



^dX 

2 dt 

dY 
dt 


= w+ {- cos X + c6^/^)Y + Z, 
+ 4y, 

T T 


dZ _ 
dt 


T T 



(9) 

where ^y and ^z are uncorrelated Wiener processes of strength 2, and 5 = Dz/Dy. Thus, 
Y and Z are completely uncorrelated, and the triple {X, Y, Z) follows a Markov process, 
with an associated FP equation. The form of Eqs. ^ is similar to models used to describe 
pendulum motion in the overdamped limit, the orientation of electric dipoles in an external 
field, the Josephson junction [T71 dB], and semiclassical lasers [13]. Note also that Eq. ^ 
contains limiting cases that depend on the value of r. For r — >■ cx3, we obtain the adiabatic 
regime described by Lapeyre [7]. Instead, the focus of this paper is on rapidly-forced regimes, 
wherein r is either comparable to the root-mean-square values of the fluctuations, or tends to 
zero; the latter case gives rise to some analytical results concerning the orientation dynamics. 
We compute the PDF of angles X, and growth rates A. These are marginal PDFs that 
can be derived from the PDF P {X, Y, Z, t), which satisfies the following FP equation: 



where 



^ou — -7^ 
Tdy 



f-^-^-/x(^^)' 


(10) 


Dy d^ 1 d Dzil-c') d^ 
r2 dY^ TdZ^^ ' r2 dZ^ 


(11) 



is the OU operator associated with the X — Y variables, and 

V = 2[w+{-cosX + c5^/2) Y + Z] 

is the drift velocity. We solve the stationary FP equation with X G [— 7r,7r] and Y,Z & 
(— oo, oo). The PDF of angles is obtained from the stationary solution P (X, F, Z): 

/»00 /"OO 

Px{X)= / dY / dZP(X,F,Z). 
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Figure 1: Marginal probability distributions obtained from the stationary FP equation, (a) The 
distribution in X-Y space; (b) The distribution in Y-Z space. The latter is Gaussian, with decay 



scales in agreement with the theoretical values in Eq. (12) 



Since X and Y are OU processes, we also have 

Pyz (F, Z) = r dXP (X, Y, Z) oc e-^'-/(2^>')e-^'^/[2^-(i-^')] . (12) 

Finally, the PDF of growth rates A = — 2y sinX is computed through a coordinate trans- 
formation of the F- variable 



^a(A) 



P 



XY 



X 



A 



' -2sinX/ 2|sinX 



1 



rdX, 



(13) 



and this transformation is legitimate because the Jacobian diverges only where P vanishes, 
and P vanishes rapidly as F — )■ ±oo. 

Sample numerical results concerning the marginal PDFs are shown in Figs. [lH2J The 



diffusion equation ( 10 ) is solved using semi-implicit spectral method; as a validation of our 



technique, we compute the marginal PDF Pyz- This is Gaussian, with decay scales in 



agreement with Eq. (12). Thus, we are satisfied with the accuracy of the method. The PDF 



of angles X is computed from the full solution P (X, F, Z) using numerical integration: this 
possesses two maxima to the right of X = ±7r/2. The shift in the maxima away from ±7r/2 
is due to the finite value of c in the calculation, while the non-invariance of the curve under 
X — > —X is due to the finite value of the drift w. When these parameters are set to zero, 
alignment with the directions ±7r/2 is achieved (Sees. IV -Vl). Finally, the mean growth rate 
of the tracer gradient is positive, since the distribution of growth rates is skewed (Fig.|2](b)). 
This suggests the favoured orientation is close to the compressive direction. This situation 
corresponds to that described by Gonzalez and Paranthoen |TDj, where in rapidly- varying 
flows, the tracer gradient aligns with a direction flxed by (r), which coincides with the 
compressive direction only when (r) = 0. 

Lastly, we focus on the limiting case where r — ?• 0. In this case, an explicit formula for 
the PDF of angles exists (up to quadratures). This enables us to prove some rigorous results 
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Figure 2: (a) The PDF of angles according to the SDK (Is]), for nonzero values of w and c. The 
favoured orientations are shifted away from =b7r/2; (b) The PDF of growth rates. A clear preference 
for positive growth is present, indicating that the tracer gradient prefers an orientation close to 
the compressive direction. Here w = c = 0.5 and Dy = Dz = 1. 



concerning the orientation of the tracer gradient relative to the compressive direction; it 
also facilitates a detailed parameter study to investigate whether the mean growth rate is 
always positive. Now in the Wiener limit, the stationary PDF of angles satisfies the ordinary 
differential equation 

J = A{X)P{X)-^ {P (X) B (X)) , (14) 

where J is the probability current and in the Stratonovich interpretation, the coefficients A 
and B have the following form: 



where g (X) 



AiX) = w + Dyg'iX)[giX) + c5'/'], 
B (X) = Dz{l- c') + Dy [g (X) + c6'/'] ' 

cosX. Introducing the effective potential 

"^ A (X') 



W(X) 



.BiX' 



-dX', 



P(X) 



-UiX) 



B{X) 



N-J e"(^')dX' 



where N and J are constants of integration fixed by normalisation and by the periodicity 
condition 

fi(-7r)P(-7r) [e"(--)-e"W] 



J 



r e"(^)dX 

J— TT 



(15) 



where w 7^ (the case w = is discussed below). Thus, P (X) is 27r-periodic on [— vr, vr 
[T7] . An expression for the mean growth rate (A) = — (2ysinX) now follows from an 




(a) w = 



^ 1 



-c=0 
-c=0.5 
• c=0.95 



(b)w==0.5 



^ 1 





i 


c=0 

- - -c=0.5 
• - c=0.95 


\ 





10" 10' 

8 



{c)'w — 5 



Figure 3: Parameter study of the mean growth rate (A): it is positive for all positive values of S. 
For w = 0, (A) approaches zero as (5 — t- 0. For non-zero if-values, (A) approaches positive limiting 
values as (5 —7- and as (5 — t- oo. 

application of the Forutsu-Novikov theorem [T71 [19]: 

{Y it) (?(") (X)) = 4Dy r dt'6 it - t') (^g^"^'^ {X {t')) [g {X {t')) + c6^''] ). 



Setting n= 1 gives and recalling that Aq = 0, we obtain the identity 

(A) = SDy ((cos^X) -C(5^/2(cosX)) . 



(16) 



Two analytical results arise from this formalism. First, when c = 0, the mean growth rate 



is definitely non-negative. Moreover, if w 

-1/2 



c = 0, it follows that J = 0. Thus, Px oc 
i?(X)~^^^, and Px attains its maximum at X = ±7r/2. In other words, for uncorrelated 
forcing with zero mean vorticity, the tracer aligns with the compressive direction, on average. 
In this case, 



(A) = 8Dy5 



1 + 






-1 



(17) 



where E and F are incomplete elliptic functions of the first and second kinds, respectively. 
To study the general case, we turn to the numerical results in Fig. |3} where we focus on 
three parameters, w/Dy, S, and c (we set Dy = 1)- The mean growth rate (A) is positive 
in all of the parameter studies. For w = 0, it tends to zero as 5 — )■ 0. For other w-values, 
it approaches definite, positive limiting values as 5 — )■ and 5 — ;■ oo. The zero-limit is c- 
independent, while the (5 — ?■ oo limit depends on c: the stronger the correlation, the smaller 
the asymptotic value taken by (A). We now compare the analytical results of the model 
with numerically-simulated flows. 



IV. THE RANDOM-PHASE SINE FLOW 

In this section we compute the orientation statistics for the random- phase sine flow [31 EH] . 
This flow is amenable to a comparison with our stochastic model because the randomisation 
of the phases breaks the invariant tori, induces homogeneity, and promotes mixing [21] . This 
is in contrast to the unrandomised sine flow studied elsewhere [TOl. 
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The flow is quasi-periodic and is defined as follows. In the j*^ period, 

u = Aq sin {ky + (pj) , v = 0, (18a) 

for the first half-period, and 

M = 0, V = Aosin^kx + ipj) (18b) 

for the second. Here r is the period of the flow, Aq is the flow amplitude, and (pj and ipj are 
random phases. 

In the j^^ cycle of the sine flow, the strain and the vorticity have the form 

d = \AakHr it) cos {ky + (pj) + \Aq (1 - Hr (t)) cos {kx + tfjj) , (19a) 

u = AokHr (t) cos {ky + (f)j) + Aok (1 - Hr (t)) cos {kx + ipj) , (19b) 

where H^ (t) is equal to unity in the first half-period and is zero in the second. Equations ( 19 ) 
imply that d and u are not correlated with each other, and that their autocorrelation func- 
tions are nonzero only in a half-period, wherein they take constant values. Thus, to mimic 



the stirring protocol (18) by an OU process, we identify the period r of the sine flow with 
the OU decay timescale, and take w = c = and Dy = Dz = A^k"^ / {2t) in Eq. (|9]). 
Furthermore, since if = Const, for the sine flow, 

{dt + u -V) X = uj- 2dcosX, (20) 



where X = 2 (^f3 — ip + ^tt) = 2/3. Using Eqs. (19), this becomes 

dX 



1 -h cos X 
in the first half-period, and 

dX 



-Aok cos {kyn + (pn) (21) 



Aok cos {kxn+i + ipn) (22) 



1 — cos X 

in the second. The growth rate of the tracer gradient along a trajectory originating at Xq is 

thus 

1 ^ 
A (a^o) = lim - V \-d''^^'^ sin X^+^Z^ _ ^n+i ^^^ ^n+ii _ /23) 

n=0 

We examine the PDF of angles X as generated by the sine flow and compare it with the 
PDF obtained from the model OU process (Fig. |4]). The results are qualitatively similar: in 
both cases, the PDF is symmetric about X = 0, with maxima close to X = ±7r/2. The sine- 
flow maxima deviate slightly from these values: this is probably due to the slow convergence 
of the PDF as the number of particles in the ensemble is increased. In the OU case, the FP 
equation gives rise to a fully-converged PDF, whose maxima occur at exactly X = ±7r/2. 
This qualitative similarity between the two PDFs is good, in view of the radically different 
noise-generating processes in each model: in the OU case, the noise terms are unbounded, 
while in the sine-flow case the fluctuations in Y and Z are bounded, |y|, |Z| < A^k. We also 
compute the PDF associated with the growth rate (Fig. Isl). This pair of PDFs for the sine 
flow demonstrates that on average, the tracer gradient aligns with the compressive direction 
of the flow. Finally, we examine the Eulerian structure of the flow by plotting the growth 
rate A (xq) (Fig. [6]). This is spatially homogeneous, thus underscoring the efficient nature of 
the random-phase sine flow in mixing the passive tracer. It also provides further justification 
for our application of the OU model to the problem in hand, since the flow statistics are, on 
average, the same along all trajectories. 
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Figure 4: PDF of the angle X as a function of flow amplitude ^o^ comparison between the sine 
flow and the OU model. The comparison yields good qualitative agreement, despite the difference 
between the statistics of the underlying forcing terms in each model. 





(a) Aa = 0.5 



(b) ^0 = 1 



(c)Ao 



Figure 5: PDF of the growth rate A as a function of flow amplitude Aq. In the main figures, we 
show a comparison between the OU model and the sine flow, while the insets contain the PDF of 
the OU model over an extended range. In the sine-flow case, the PDF takes values in a bounded 
interval, while the tails of the OU model extend to ±00. 






(a) Ao = 0.5 



(b)^o = l 
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Figure 6: The growth rate of the tracer gradient along trajectories, according to Eq. (23) 
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Figure 7: Orientation dynamics of forced two-dimensional turbulence in a statistically steady state: 
shapshots of (a) the angle (/?; (b) the angle X; (c) the strain rate d. 

V. FORCED TWO-DIMENSIONAL TURBULENCE 



In this section we compare the OU model of the orientation dynamics with the results 
of numerical simulations of forced two-dimensional turbulence. We present results for the 
solution of the vorticity equation 






+ u ■ Vw 



-1)^ VpS/^U + Q — UqUJ, 



(24) 



where Q is a forcing term and —i>ou is a damping term that prevents a buildup of energy 
at large scales. (Simulations for Rayleigh damping [22j have also been carried out, and 
yield similar results.) To obtain a stochastic driving force on a particular scale kg, we 



use the method developed by Lilly [23] (see also [IS]). The combination of damping and 
forcing in Eq. ( [24) ) yields a statistically-steady state. Moreover, the tracer, driven by the 
flow (24), also reaches a steady state, in the sense that 6'/||6'||2 exhibits the so-called strange 
eigenmode [31 El] . Thus, it is legitimate to regard the late-time statistics of the angle (3 as 
being drawn from a stationary distribution. 

To map the flow problem onto the OU dynamics, we examine the Eulerian fields 



w 



\u + dti^ + u ■ Vv?, A = sign {d) \/d? + s^. 



A snapshot of these fields at the same time is shown in Fig. [T) We solve Eq. (24) using 
a standard, semi-implicit pseudospectral method. We take the order of the viscosity to be 
p = 8; we also take z/p = 5.9 x 10"'^° in both the u- and the ^-equations, and z/q = 0.05. 
The spatial resolution is A^^ = 256^ and the timestep is set to At = 10^^. We identify 
the decay time r with Uq^ and examine Eulerian time series of w and A over the steady- 
state period of the numerical integration. These scalar variables fluctuate around certain 
mean values, such that (||A||2)t = 0.05/r, (||w||2)i = 0.94/r, and (/d^xAw)j = 0.00 (the 
brackets {■)t denote temporal averaging over a long steady-state interval). Thus, we model 
the orientation dynamics as OU processes with parameters Dy = 0.05/r, Dz = 0.94/r, 
c = 0, and w = 0. Making these identifications, we compare the PDF of the angle X and 
the growth rate A, as generated both by the model OU process, and the two-dimensional 
turbulence. The results are shown in Fig. [8] The angle PDFs for the turbulence and for 
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Figure 8: (a) The PDF of the angle X according to the 2d turbulence simulation (solid line), and 
the OU model (broken line) (b) The PDF of the growth rate. The model is in excellent qualitative 
agreement with the turbulence. 



the OU model are in very close agreement, with large maxima at exactly X = ±7r/2. The 
PDFs for the growth rate are in close qualitative agreement: the curve is asymmetric around 
A = 0, and is skewed towards positive values. Both curves fall off sharply away from A = 0. 
The first moment of the turbulence PDF is less than that of the OU process [Ai^ou = 0.15; 
Ai,TURB = 0.08), and the tail of the turbulence PDF is less 'fat' than its OU counterpart. 
Nevertheless, the close qualitative agreement between the model and the flow both in this 



section and in Sec. |IV| confirms the usefulness of the OU description in modelling the growth 
rate of a tracer gradient under rapid stirring. 



VI. CONCLUSIONS 



We have formulated a stochastic model of the orientation dynamics of a tracer gradient in 
a manner designed to mimic mixing under externally-imposed rapidly- varying stirring in two 
dimensions. The model consists of a set of SDEs for the angle (, and for the forces that act 
thereon. The model is a general one that includes the adiabatic description of Lapeyre [7j 
as a special case. However, we focus entirely on the rapid case: we analyse the SDEs using 
the Fokker-Planck (FP) equation and compute the distrubtion of orientation angles and 
growth rates. The use of the FP equation gives rise to some analytical results in the limit 
of zero correlation time, and in the special case when the mean vorticity and the correlation 
between the forcing terms both vanish: the tracer gradient aligns with the compressional 
direction, and the mean growth rate of the tracer gradient is positive and given by an explicit 
expression in 6 = Dz/Dy involving elliptic functions. For finite correlation times, numerical 
solutions suggest that the PDF is always skewed to positive values, giving rise to a positive 
mean growth rate. We compare our model with two rapid-fiow protocols: the random-phase 
sine fiow, and forced two-dimensional turbulence. The qualitative agreement between the 
model PDFs and the PDFs generated from the fiow protocols is excellent, and confirms 
the validity of the approach. Since our stochastic model is amenable both to mathematical 
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and numerical analysis (the latter in minutes on a desktop computer), we anticipate its 
generalisation to three-dimensional flows. 
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